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There has been considerable interest recently in a variety of liquids composed of line- 
like objects aligned on average along a common direction. Flux lines in high-temperature 
superconductors, strings of electric dipoles in electrorheological fluids and polymer nemat- 
ics are all systems with this common morphology. Despite the vast literature on nematic 
liquid crystals, the number of theoretical and experimental studies of polymer nematics 
is relatively small. Building on work by Meyer and de Gennes [H, statistical me- 
chanical treatments of polymer nematics have predicted, in particular, the structure 
function in the semi-dilute regime. Ao, Wen and Meyer [|5| have studied the structure of 
poly-7-benzyl-glutamate in the nematic phase by X-ray scattering. However because the 
experiments cannot reach very small momentum q, due to the forward scattering beam 
stop, they cannot test many of the main theoretical predictions which are good at these 
small q values. 

We have simulated polymer nematics with the hope of gaining some additional infor- 
mation in the behavior of the polymer nematics. In our simulation nematic order arises 
only due to the polymer stiffness - there is no microscopic, steric nematic interaction. 
We are particularly interested in the small q region of the static structure function. The 
resulting equilibrium polymer configurations could be used to compare to actual polymer 
configurations in nematic solutions of actin 0. We contrast this study to recent work 
which studied the phase behavior of 10-mers made of rigid nematogens. In our simulation 
nematic order arises only due to the polymer stiffness - there is no microscopic, steric 
nematic interaction. Because of the long length scales and slow relaxation times which 
are inherent in any polymer problem, it is simply not feasible to model in great detail the 
chemical complexity of any of the experimental polymer nematics such as actin. Fortu- 
nately since we are interested in understanding some of the basic properties of polymer 
nematics, particularly at small q (large distances), it is not necessary to include the chem- 
ical detail explicitly. For this reason, we chose to study a coarse grained pearl-necklace 
model similar to that previously studied for flexible polymer systems [§. The stiffness of 
the chain was controlled by adding a three-body bending term which could be tuned to 
vary the persistence length of our polymer nematic relative to its chain length. Using this 
model, we are able to test several of the theoretical predictions of . 

We studied an ensemble of M polymer nematics in a periodic rectangular cell. Each 
monomer of the polymer was modeled by a hard sphere of diameter a connected to its 
neighbors along the chain with a fixed bond of length a = 1.25cr. Results presented here 
were for M = 200 chains of length = 50. Since we are only interested in the static 
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equilibrium properties of the system and not their dynamics, we used a shthering snake 



algorithm which is very efficient particularly when the persistence length is large ||1C1|| . 
The algorithm models the reptation motion of a snake by attempting to remove the head 
(tail) of the chain and placing it at the tail (head) with a random bond angle. We modeled 
the bond energy as —kcos9 where 9 is the bond angle, and k is the bending stiffness 
measured in units of kj^T. If the move violates the excluded volume constraint on the 
beads, the move is rejected. If not, then the Monte Carlo move is accepted according to 
the relative energies of the original and new bond angle. The results presented here are 
from runs made using a force-biased algorithm [l^] in which moves are only rejected due to 
excluded volume. While more moves were accepted compared to a standard Monte Carlo 
algorithm, the additional complexity did not lead to significant time reduction at high 
density on either the Silicon Graphics Power Challenge or the IBM Power PC604 (runs 
were done on both). Since we were interested in rather stiff chains and the rejection rate for 
the moves was quite high, we quote all times in terms of successful moves [|12[ . Because of 
the long runs we concentrated on three areal number densities pocr^ = 0.142 (areal fraction, 
po7Ta'^/4 = 0.11), 0.32 (0.25), and 0.569 (0.45). At the highest density the system did not 
equilibrate and presumably became glassy. We will not report here on the highest density 
run. We found that the autocorrelation of the end-to-end distance for a single polymer at 
the bending stiffness we are interested in (k = 40) decayed exponentially with a relaxation 
time of about 1.5 x 10^ Monte Carlo moves. Thus we believe that after 1.2 x 10^ steps 
the two lower density runs should have equilibrated. The simulations were started from 
an initial condition in which all 200 polymers chains were aligned along the 2;-axis but 
randomly placed in the xy plane at random heights along the z-axis. At poc'^ = 0.142, 
each monomer, on average, moved 15a laterally, while at poc^ = 0.32 each moved lOcr. 
The chains were placed in a periodic rectangular box with a fixed height = Na = 62.5cr 
in the 2;-direction and a variable xy area. For the densities we ran at, we found a Maier- 
Saupe order parameter {S = [3(cos2^) - l]/2) of 0.049 and 0.849 for poa^ = 0.142 and 
0.32, respectively, where 9 is the angle between the polymer bond and the principal axis 
with the largest eigenvalue of the bond-moment tensor Tij = (titj), where t is the vector 
pointing along a bond. Since there was no explicit nematic field added, there is nothing to 
stop the ordered domain from rotating to a new direction. In the denser run we found that 
the nematic axis differed from the f-axis by approximately 7°. In the most dilute system, 
since there was no nematic order, the issue was moot. We thus take the z-direction to be 
the ordering direction for the analysis. In addition, since the starting states had S = 1, 



2 



the lower average values of the order parameter give us an independent confirmation that 
the system has equilibrated. 

We calibrated the bending energy k with polymer persistence length for a dilute chain. 
The persistence length was determined from the bond angle-bond angle correlation length. 
In Fig. 1, we show the persistence length Lp/a as a function of k calculated from the 
first half (solid line) and second half (dashed line) of a 2.4 x 10^ step run. Because of 
the excluded volume interaction, a monomer cannot bend back on its neighbor and the 
maximum bending angle of a single bond is 133°. This is why the persistence length Lp/a 
is greater than 1 for k = 0. In addition, note that if we were to take the continuum limit, 
the bending energy would become |k J{ds9)^- This simple theory leads to a persistence 
length Lp/a = n/k^T, in agreement with our calibration when the stiffness is such that 
self-avoidance is not an issue. In an attempt to match to the experiments of Ao, Wen and 
Meyer ||^ where the polymer persistence length was 80% of the chain length and what can 
be run in a reasonable amount of time, we concentrate on k = 40 for the rest of this paper. 
This corresponds to a persistence length Lp/a ~ 39 and Lp/[{N — l)a] = 80%. 

Let us recall the theoretical expectations for the structure function 

S{q±,qz) = 2 ' (1) 

Po 

where q± is the magnitude of q± = {qx: qy^O), P = Po + Sp is the local areal number density 
and Po is the average mean density. For infinitely long polymers, we expect no scattering 
along the qz axis due to the conservation of polymer density along the z-axis. As described 
by Taratura and Meyer [|13[ , the projection of the polymer tangent into the xy-plane t and 



the areal density p satisfy 

dzP + V±-pt = (2) 

When (jj^ = (H) implies that dzp{q± = O^z) = 0. Hence there is no density contrast and 
there will be no scattering along the q± = axis. When the polymers are finite in extent, 
(^ is modified by adding sources and sinks for polymer heads and tails to the right hand 
side of the equation. This will lead to some scattering along the qz axis Q controlled 
by the typical polymer length. In addition we expect Bragg-like peaks on the q± axis, 
corresponding to the incipient columnar crystal being formed by the polymers. Together 
these two predictions suggest the typical bowtie anisotropic pattern of X-ray scattering. 



The average static scattering function 5'(q) can easily be determined for the simulated 
system from 



We took the average every 1.2 x 10 moves. Because of the periodic boundary conditions 
each Qi has to be commensurate with the dimensions of the cell. This means that the 
smallest nonzero Qi = 2n/Li, where i is one of the Cartesian coordinates. Thus as po 
increases at fixed M, the smallest accessible qi increases, which is one reason we did not 
study higher density systems. In Fig. 2, we present the results for the two dimensional 
structure factor S{q±, qz) for the two densities studied. 

Our data agree qualitatively with the first expectation and quantitatively with the 
second. In addition we calculated the average polymer persistence length for the densities 
considered. We find Lp/a = 39 {poa'^ = 0.142) and 95 (0.32). This is consistent with the 
increased amount of nematic ordering with increasing density. Note that the value at the 
lowest density is close to the dilute value of 39. 

While the structure function data in the {q±, qz) plane does not have enough resolution 
for the fitting of the full two dimensional surface, some information may be gleaned from 
collapsing the data onto the q± axis. Indeed, the derived function 



is the structure function of the polymer nematic in any fixed z cross section 0] . Since the 
numerically determined structure function is only computed for Qz = — IOOtt/L^ . . . lOOn/Lz 
(where Lz = Na = 62.5a is the height of the box) the sum performed to calculate S2{q±) 
will be somewhat smaller than the actual value of S2{q±)- We see the formation of true 
liquid-like structure in the nematically ordered system and of almost gas-like order in the 
isotropic system. This is certainly reasonable: the isotropic system is in no sense a directed 
line-liquid and we should not expect any correlated behavior in a constant 2;-slice. Note 
that in the ordered run S2{q±) has maxima at multiples of go = 47r/(v^ao) where gq is 
the average inter-polymer spacing, indicative of the incipient crystalline order. 

In this study we have presented the first Monte Carlo study of a nematic polymer 
liquid crystal for densities in the semi-dilute regime. We found that for intermediate 
densities, the static structure function is anisotropic and has the predicted bowtie shape 
in the (gj^, g^) plane. Using a sum-rule, similar to that relating the structure function at 
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zero momentum to the bulk compressibility ||^] we could extend this work to study the 
scaling of the bulk modulus over a range of densities. In addition, the sum rule could be 
checked differently by varying the average polymer length. Polydispersity should not effect 
the theoretical predictions: this could be checked as well. Finally, a study which varied 
the persistence length at fixed density would also be enlightening and could be compared 
to theory. 
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Figure Captions 

Fig. 1. Single polymer persistence length (in units of the center to center spacing 
o = 1.25a) as a function of the stiffness parameter k. The solid line shows 
the persistence lengths measured from a 1.2 x 10^ step run, while the dashed line 
shows the persistence lengths measured in the second half of a 2.4 x 10^ step run. 

Fig. 2. Density plots of the polymer structure function in the {q±, q^) plane. Here go = 
47r/-\/3aO) where oq is the average inter polymer spacing. Note that there is little 
scattering along q± = 0, in agreement with theoretical arguments, (a) areal 
number density 0.142 (b) areal number density 0.32. 

Fig. 3. Two-dimensional in-plane structure function for areal number densities 0.142 and 
0.32. This is obtained by S2{q±) = J ^S{q±,qz). The curves are successively 
shifted along the intensity axis for clarity. 
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